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Properties of the Olami-Feder-Christensen (OFC) model of earthquakes are studied by numerical 
simulations. The previous study indicated that the model exhibited "asperity" -like phenomena, i.e., 
the same region ruptures many times near periodically [T. Kotani et al, Phys. Rev. E 77, 010102 
(2008)]. Such periodic or characteristic features apparently coexist with power- law- like critical 
features, e.g., the Gutenberg- Richter law observed in the size distribution. In order to clarify the 
origin and the nature of the asperity-like phenomena, we investigate here the properties of the 
OFC model with emphasis on its stress distribution. It is found that the asperity formation is 
accompanied by self-organization of the highly concentrated stress state. Such stress organization 
naturally provides the mechanism underlying our observation that a series of asperity events repeat 
with a common epicenter site and with a common period solely determined by the transmission 
parameter of the model. Asperity events tend to cluster both in time and in space. 



INTRODUCTION 

Statistical properties of earthquakes have attracted 
much interest in seismology as well as in statistical 
physics pj, |2]. Statistical properties of earthquakes are 
often characterized by power-laws, e.g.^ the Gutenberg- 
Richter (GR) law observed in earthquake size distribu- 
tion, or the Omori law observed in the time evolution 
of aftershocks frequency. The concept of "self-organized 
criticality" (SOC) was advocated by P. Bak In this 
view, the Earth's crust is supposed to be in a critical 
state, while seismicity has a close analogy to critical phe- 
nomena of thermodynamic second-order transition. 

Several statistical models of earthquakes embodying 
such SOC features of earthquakes have been proposed 
and studied. One standard model is the so-called spring- 
block or the Burridge-Knopoff (BK) model, in which an 
earthquake fault is modeled as an assembly of blocks 
mutually connected via elastic springs which are slowly 
driven by external force [4-6]. While this model some- 
times exhibits critical or near-critical properties under 
restricted conditions, generic properties of this model 
are not necessarily critical, often exhibiting "character- 
istic" properties accompanied by characteristic energy 
and time scales. For example, the BK model with the 
nearest-neighbor 0, HI or the long-range interaction ^ 
has turned out to exhibit either "supercritical" , "subcrit- 
ical" or "near-critical" behavior depending on the param- 
eters of the model. 

Another statistical model extensively studied in sta- 
tistical physics in the context of SOC might be the so- 
called OFC model, which was first introduced by Olami, 
Feder and Christensen (OFC) as a further simplification 
of the BK model [10]. It is a two-dimensional lattice 
model where the rupture propagates from lattice site to 
its nearest-neighbor sites in a non-conservative manner, 
often causing multisite seismic events or "avalanches". 
Numerical studies have revealed that the OFC model 
exhibits apparently critical properties such as the GR 



law lol-lll or the Omori law [HI, although there still 
remains controversy concerning whether the behavior of 
this model is strictly critical jl2| or only approximately so 
[l3-[l6|. In any case, the OFC model has been regarded 
as a typical non-conservative model exhibiting SOC. 

Earthquakes in nature sometimes exhibit characteris- 
tic properties accompanied by characteristic energy and 
time scales [T, 'i]. In contrast to its critical feature, pos- 
sible "characteristic" feature of the OFC model, if any, 
has attracted much less attention so far. Botani and 
Delamote observed that synchronized clusters of vari- 
ous sizes existed in the steady state of the OFC model 
with a large, but finite lifetime [l^- More recently, 
Kotani, Yoshino and Kawamura have revealed strikingly 
characteristic features of the model, by demonstrating 
that the local recurrence-time distribution of the OFC 
model exhibits a sharp 5-function-like peak correspond- 
ing to rhythmic recurrence of events with a fixed period 
uniquely determined by the transmission parameter of 
the model [18], together with a power-law-like tail corre- 
sponding to scale-free recurrence of events. In fact, the 
OFC model was found to exhibit phenomena closely re- 
sembling the "asperity" known in seismology. Such an 
asperity repeats to rupture sequentially, often more than 
ten times, with almost the same period and with the 
same site as an epicenter (triggering site) [18]. Statistics 
of epicenter sequences of the model was also studied by 
Peixoto and Prado in the framework of a growing com- 
plex network [Toj . 

In the present paper, following Ref.fis'], we further in- 
vestigate spatiotemporal correlations of the OFC model, 
with particular emphasis on its characteristic features. 
In particular, we wish to unravel more detailed proper- 
ties of the "asperity" involved in the model. It is found 
that, in the formation of the asperity, self-organization 
process occurs in which the stress in the asperity region 
gets more and more "quantized" to discretized values. 
This self-organization process promotes an epicenter of 
the preceding event to become an epicenter of the next 
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event again with the period uniquely determined by the 
transmission parameter of the model. It is also observed 
that the epicenter of large events tends to lie at the tip 
(or at the corner) of the rupture zone, rather than in its 
interior. 

After defining the model and explaining some of its 
basic properties in §2, we present our numerical data of 
the local recurrence-time distribution and of the stress 
distribution in §3. Various asperity characteristics are 
investigated in detail in §4. The relation between the 
self-organization of the highly concentrated stress state 
and the asperity characteristics are clarified. Section 5 
is devoted to summary and discussion. Analysis of the 
time evolution of the stress concentration is given in Ap- 
pendix. 

THE MODEL 

In the OFC model, "stress" variable fi {fi > 0) is as- 
signed to each site on a square lattice with L x L sites. 
Initially, a random value in the interval [0,1] is assigned 
to each /i, while fi is increased with a constant rate uni- 
formly over the lattice until, at a certain site z, the fi 
value reaches a threshold, /c = 1. Then, the site i "top- 
ples", and a fraction of stress afi {0 < a < 0.25) is 
transmitted to its four nearest neighbors, while fi itself 
is reset to zero. If the stress of some of the neighboring 
sites j exceeds the threshold, i.e., fj > /c = 1, the site 
j also topples, distributing a fraction of stress afj to its 
four nearest neighbors. Such a sequence of toppHngs con- 
tinues until the stress of all sites on the lattice becomes 
smaller than the threshold /c. A sequence of toppling 
events, which is assumed to occur instantaneously, cor- 
responds to one seismic event or an "avalanche". After 
an avalanche, the system goes into an interseismic period 
where uniform loading of / is resumed, until some of the 
sites reach the threshold and the next avalanche starts. 

The transmission parameter a measures the extent of 
non-conservation of the model. The system is conserva- 
tive for a = 0.25, and is non-conservative for a < 0.25. 
A unit of time t is taken to be the time required to load 
/ from zero to unity. 

In the OFC model, boundary conditions play a cru- 
cial role. For example, SOC state is realized under open 
boundary conditions, but is not realized under periodic 
boundary conditions. The model under open boundary 
conditions goes into a special transient state where events 
of size 1 (single-site events) repeat periodically with pe- 
riod 1 — 4a These single-site events occur in turn 
in a spatially random manner, but after time 1 — 4a, 
the same site topples repeatedly. Although such a pe- 
riodic state consisting of single-site events is a steady 
state under periodic boundary conditions, it is not a 
steady state under open boundary conditions because 
of the boundary. Indeed, clusters are formed near the 



boundary, within which the stress values are more or less 
uniform, and gradually invades the interior destroying 
the periodic state. Eventually, such clusters span the en- 
tire lattice, giving rise to an SOC-hke steady state. Such 
clusters might be formed via synchronization between the 
interior sites and the boundary sites, the latter having a 
slower effective loading rate due to the boundary [20]. 

In our simulation, the lattice studied contains 256 x 
256 sites with open boundary conditions, the pseudo- 
sequential updating proposed by Pinho et al being uti- 
lized [21] . Total of 2 X 10^ avalanches are generated where 
initial 10^ events are discarded as transients to ensure 
that the measurement is done in the steady state. 

It has been realized that the OFC model exhibits an 
SOC-hke critical or near-critical property. This is most 
clearly exemplified in a power-law-like size distribution 
of avalanches. As an example, we show in Fig.l the size 
distribution of avalanches on a log-log plot for several 
values of the transmission parameter a. The avalanche 
size s is defined by the total number of "topples" in a 
given avalanche, which could be larger than the num- 
ber of toppled sites because multi-toppling is possible in 
a given avalanche. In fact, it is observed that multi- 
toppling rarely occurs in the model except in the conser- 
vation limit or in the regime very close to it. 
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FIG. 1: (Color online) The size distribution of seismic events 
(avalanches) for various values of the transmission parameter 
a. The slope of the data gives the value of 1 + which is 
shown in the figure. 

As can be seen from Fig.l, a straight-line behavior cor- 
responding a power-law is observed consistently with the 
previous works. The slope representing the 5- value is 
not universal varying from 0.90 to 0.22 as a is varied 
from 0.17 to 0.245. The power-law feature is weakened 
as one approaches the conservation limit. The 5- value 
is known to come around 2/3 in real seismicity, while 
the 5- value of the OFC model is distributed around this 
value [10,]. 
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Hergarten et al observed that the OFC model also 
exhibited another well-known power-law feature of seis- 
micity, i.e., the Omori law (or the inverse Omori law) 
describing the time evolution of the frequency of after- 
shocks (foreshocks) [l3|. We show in Fig. 2(a) on a log- 
log plot the frequency of aftershocks as a function of the 
time elapsed after the mainshock t. As a mainshock, we 
consider an event of its size greater than s > Sc = 100. 
Again, a straight-line behavior corresponding a power- 
law is observed consistently with Ref.[l3]. Note that Her- 
garten et al imposed free boundary conditions, while we 
impose open boundary conditions here. The slope rep- 
resenting the Omori exponent p is again not universal 
depending on the parameter a as p = 0.84, 0.69 and 0.03 
for a = 0.17,0.20 and 0.23, respectively. Since the p- value 
is known to come around unity in real seismicity, the p 
value of the OFC model is not necessarily close to real 
observation. 

Similar results are obtained also for foreshocks: See 
Fig. 2(b). Here the p- value describing the inverse Omori 
law is equal to p = 0.71, 0.50 and 0.01 for a = 0.17,0.20 
and 0.23, respectively. As one approaches the conserva- 
tion limit a = 0.25, both aftershocks and foreshocks tend 
to go away. Our results are consistent with the results 
reported by Hergarten et al [13]. 

As such, the OFC model certainly possesses critical or 
near-critical properties like the GR law and the Omori 
law. By contrast, its characteristic properties, most no- 
tably revealed in the existence of "asperity" [Tsj, have 
received much less attention. In the next section, we 
investigate such characteristic properties of the OFC 
model, which apparently coexist with the critical prop- 
erties mentioned in this section. 
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SIMULATION RESULTS ON THE 
CHARACTERISTIC PROPERTIES 

In this section, we report on our simulation results on 
the local recurrence-time distribution and the stress dis- 
tribution, which turn out to exhibit pronounced charac- 
teristic properties. 

Let us begin with the definition of the local recurrence 
time. In case of the globally defined recurrence time, the 
next avalanche to measure recurrence may occur any- 
where on the entire lattice. In view of the ordinary sense 
of earthquake recurrence, however, it might be more nat- 
ural to introduce the recurrence time T and its distri- 
bution function P{T) locally. Thus, one may define the 
local recurrence-time as the time passed until the next 
avalanche occurs with its epicenter lying in a vicinity of 
the preceding avalanche, say, within distance r (in units 
of lattice spacing) of the epicenter of the preceding event, 
rather than in some remote place far away from the epi- 
center of the preceding event. One can also introduce the 
size threshold Sc to look at the recurrence of large events 



FIG. 2: (Color online) The time dependence of the frequency 
of aftershocks (a), and of foreshocks (b), on a log- log plot for 
several values of the transmission parameter a. Mainshocks 
are the events of their size greater than s > Sc = 100. The 
time t is measured with the occurrence of a mainshock as an 
origin. Aftershocks and foreshocks are defined as events of 
arbitrary sizes which occur in the vicinity of mainshock with 
the range parameter rc = 10. 



of their size greater than s > Sc. 

The local recurrence-time distribution of the model, 
P{T)^ was calculated in Ref.[18] with varying the range 
parameter Tc. The computed P{T) exhibited a sharp S- 
function-like peak at T = T* = 1 — 4a which grew as rc 
got smaller, indicating that many (though not all) events 
of the OFC model repeated with a fixed time-interval 
T = T*. The peak position turned out to be independent 
of the range parameter Tc, the size threshold 5c, and the 
lattice size (as long as it was not too small). 

In Fig. 3, we show the local recurrence-time distri- 
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but ion P{T) for avalanches whose size is greater than 
5 > Sc = 100 for fixed Tc = 10, with varying the 
transmission parameter a toward the conservation hmit 
a = 0.25. As a is increased toward a = 0.25, the S- 
function peak is graduahy suppressed with keeping its 
position strictly at T = 1 — 4a. For a = 0.245, the S- 
function peak is no longer appreciable at the expected 
position T = 1 - 4 X 0.245=0.02: See the arrow in the 
figure. Clearly, the (5- function peak of P{T) goes away 
toward the conservation limit a = 0.25. 
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FIG. 3: (Color online) Log-log plots of the local recurrence- 
time distributions of large avalanches of their size s > Sc = 
100 for a fixed range parameter rc = 10, with varying the 
transmission parameter a. The arrow in the figure represents 
the expected peak position for a = 0.245 corresponding to 
the period T* = 1 — 4a = 0.02. The inset is a magnified view 
of the main peak for the case of a = 0.17. Log- log plots of 
the local recurrence-time distributions of large avalanches of 
their size s > Sc = 100 for a fixed range parameter rc = 10, 
with varying the transmission parameter a. The arrow in the 
figure represents the expected peak position for a = 0.245 
corresponding to the period T* — 1 — Aa — 0.02. The inset is 
a magnified view of the main peak for the case of a = 0.17. 



As can be seen from Fig. 3, in the longer time regime 
T > T^^ P{T) exhibits behaviors close to power laws 
whereas the deviation from the power law becomes 
appreciable when one approaches the conservation limit. 
Hence, in earthquake recurrence of the model, the char- 
acteristic or periodic feature, i.e., a sharp peak in P{T) 
at T = T*, and the critical feature, i.e., power-law-like 
behaviors of P(T) at T > T*, coexist. Note that, while 
the peak at T = T* is sharp, it is not infinitely sharp with 
a finite intrinsic width as shown in the inset of Fig. 3. The 
peak has an asymmetric shape with a tail on shorter-T 
side. 

Furthermore, as was already noticed in Ref.[l8|, the 
periodic events contributing to a sharp peak of P{T) 
( "peak events" ) possess a power- law- like size distribution 
very much similar to those observed for other aperiodic 
events. This is demonstrated in Fig. 4, where the size 



distribution of the peak events is given for the case of 
a = 0.23, in comparison with the one for general events 
of 5 > Sc = 100. (At this a- value, a = 0.23, a slight 
excess is discernible at larger s as was already noticed in 
Fig.l above, both for all events and for the peak events.) 
In that sense, the peak events certainly appear to be char- 
acteristic in time, but not necessarily so in their size. 
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FIG. 4: (Color online) The size distribution of the peak 
events as compared with the one of general events of s > 
Sc = 100. The transmission parameter is a = 0.23. The 
peak event is defined here as an event of its size greater than 
s > Sc = 100 with its recurrence time belonging to the main 
peak of the local recurrence-time distribution function P(T) 
given in Fig. 3. 

The fact that the peak events possess a power-law like 
broad distribution in their size s also explains our ob- 
servation of Fig. 3 that the computed P{T) exhibits a 
multi-peak structure consisting of sub-peaks located at 
multiples of T* = 1 — 4a. Reflecting such a size distribu- 
tion of peaks events, many peak events occur with their 
size above our size threshold Sc slightly. In such a situa- 
tion, the "next" event occurring at T = T* after the flrst 
one might occur with its size slightly below the thresh- 
old 5c, and the "second- next" event occurring T = 2T* 
after the flrst one might be counted as the next event in 
measuring the local recurrence time (remember that the 
peak in P{T) has a small but flnite width). In Ref.^lH, 
the major cause of the multi-peak structure of P{T) was 
ascribed to the one associated with the range parame- 
ter Tc, but we now have found that the major cause of 
the multi-peak structure of P{T) is the one associated 
with the size threshold Sc rather than with the range 
parameter Vc'. Refer also to our analysis of the asperity 
characteristics in §4. 

We have examined whether there exists appreciable 
difference between the peak events and the other events 
with regard to their aftershock/foreshock properties. 
However, both types of events accompany very much 
similar Omori-law-type aftershock/foreshock sequences. 
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Concerning their aftershock/foreshock sequences, we 
cannot identify any appreciable difference between the 
peak events and the other events. 

In Fig. 5, we show for several a- values the time- 
averaged stress distribution in the steady state D{f) av- 
eraged over all sites of the lattice. Though the system 
is loaded with a constant rate, the stress distribution is 
not uniform in the interval [0, 1]. Rather, for the case 
of a < 0.2, D{f) exhibits distinct steps. For a < 0.2, 
these steps appear at appropriate multiples of the trans- 
mission parameter a, i.e., at / = na for / < 1/2 and at 
1 — na for / > 1/2, with n being an integer. For larger 
a > 0.2, no step seems to appear at multiples of a, but 
D{f) still exhibits a noticeable structure. In this sense, 
the stress distribution also exhibits a characteristic fea- 
ture. Such a structure in D{f) is likely to be originated 
from large-scale avalanches. 
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FIG. 5: (Color online) The time- averaged stress distribution 
in the steady state averaged over all sites of the lattice for 
several values of the transmission parameter a. 



SIMULATION RESULTS ON THE ASPERITY 
CHARACTERISTICS 

In fact, the peak events, or near-periodic events, are 
closely related to the asperity-like phenomena [18]. The 
asperity represents a cluster of sites which ruptures si- 
multaneously in a given avalanche, repeating many times 
with the same epicenter site and with almost the same 
period T = T* = 1 — 4a corresponding to the peak of the 
local recurrence-time distribution function P{T). 

We show in Fig. 6 typical snapshots of the stress dis- 
tribution immediately before and after an asperity event 
for the case of a = 0.2. Within the asperity region, uni- 
form stress state is formed before the rupture giving rise 
to a synchronized simultaneous rupture. Then, a discon- 
tinuous drop of the stress within a rupture zone (asper- 
ity) occurs before and after the avalanche. In fact, the 
same cluster except for a minor difference ruptures again 



after the time T*. Such a rhythmic recurrence of rup- 
ture often repeats more than ten times. There is a clear 
tendency that such "repeating times" of the asperity se- 
quence tend to be greater for smaller a. In particular, the 
asperity-like phenomena themselves hardly exist near the 
conservation limit (remember that the local recurrence- 
time distribution P{T) shown in Fig. 3 no longer exhibits 
a peak at a = 0.245). 




FIG. 6: (Color online) Snapshots of the stress distribution 
for the case of a = 0.2; (a) immediately before a large event 
at time t = to, (b) immediately after this event, (c) im- 
mediately before the following event which occurs at time 
t = to + T*(T* = 0.2), and (d) immediately after this second 
event. Two events are of size s — 15891 and s — 15910 on 
a L = 256 lattice. The region surrounded by red bold lines 
represents the rupture zone, while the star symbol represents 
an epicenter site which is located at the tip of the rupture 
zone. 

The stress distribution in the asperity region tends to 
be "discretized" to certain values. In Fig. 7(a) and (b), 
we show for the case of a = 0.17 the stress distribution 
D{f) of the asperity sites immediately before (a) and af- 
ter (b) an avalanche, averaged over asperity events. We 
define here the asperity as a ruptured cluster of its size 
s > Sc = 100 belonging to the peak event of the local 
recurrence-time distribution function. As can be seen 
from the figure, D{f) now consists of several "spikes" lo- 
cated at appropriate multiples of the transmission param- 
eter Of, i.e., at 1 — na before the rupture and at / = na 
after the rupture, with n being an integer. One can now 
see that the structure found in the time-averaged stress 
distribution function shown in Fig. 5 is originated from 
those spikes associated with the asperity event. 

In fact, such a tendency of "quantization" or "spikes" 
also exists in the non-asperity events as well, but in a 
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FIG. 7: (Color online) The stress distribution D{f) of each 
site contained in the rupture zone of the asperity event, just 
before (a) and after (b) the asperity event. An asperity event 
is defined here as an event of its size greater than s > Sc — 
100 belonging to the main peak of the local recurrence-time 
distribution function. The transmission parameter is a = 
0.17. The inset is a magnified view of the main peak. 



less pronounced manner. In Fig. 8, we show the stress 
distribution at the time of toppling of each site contained 
in the rupture zone for both cases of the asperity events 
and the non- asperity events. The stress at the time of 
toppling generally exceeds the threshold value /c = 1 
(except for an epicenter site where it is identically unity) . 
As can be seen from Fig. 8, the stress value at the time of 
toppling is more concentrated on the threshold value fc = 
1 for the asperity events than that for the non-asperity 
events, the latter exhibiting a broader tail toward larger 
values of / > 1. 

Furthermore, as the asperity events repeat, the ten- 
dency of the stress concentration is more and more en- 
hanced. In Fig. 9, we show the time sequence of the 
stress distribution at the time of toppling for the asperity 
events. As the asperity events repeat, the stress distri- 
bution tends to be narrower, being more concentrated on 
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FIG. 8: (Color online) The stress distribution D{f) at the 
time of toppling of each site contained in the rupture zone of 
large events with s > Sc = 100. The data are given for both 
cases of the asperity events and the non-asperity events. An 
asperity event is defined here as an event of its size greater 
than s > Sc = 100 belonging to the main peak of the lo- 
cal recurrence-time distribution function. The transmission 
parameter is a = 0.17. 



the threshold value /c = 1. In fact, one can show that 
the stress distribution at the time of toppling tends to be 
more concentrated on the threshold value /c = 1 as the 
asperity events repeat. Namely, once each site starts to 
topple more or less at similar stress values close to the 
threshold value /c = 1, this tendency is more and more 
evolved as the asperity events repeat. The stress concen- 
tration tends to be self-organized. In such a situation, 
all sites in an asperity cluster (except for an epicenter 
site) possess the stress value close to 1 — a just before 
the event starts. Assume that an arbitrary site involved 
in the interior of the asperity cluster possesses the stress 
value f = 1 — a -\- X with x << 1, just before the event 
starts. One can show that, when this event is over and 
the interseismic period of T* = 1 — 4a has passed, the 
stress value at this site becomes l — a-\-S*x with S* being 
less than unity under rather general conditions, as long 
as the stress values at neighboring sites just before the 
event are also close to 1 — a. Details of the derivation are 
given in Appendix. 

We note that, if all sites topple at the stress value 
close to the threshold /c = 1 in the asperity events, it 
immediately explains our observation that the period of 
the asperity events, corresponding to the main peak of 
the local recurrence-time distribution, is equal to T = 
T* = 1 — 4a. To see this, one only needs to remember the 
conservation law of the stress, z.e., the stress dissipated 
at the time of toppling, which is 1 — 4a per site if the 
toppling occurs exactly at / = 1, should match the stress 
loaded during the interval time T. 
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FIG. 9: (Color online) The time sequence of the stress distri- 
bution D(f) at the time of toppling of each site contained in 
the rupture zone of the asperity events. An asperity event is 
defined here as an event of its size greater than s > Sc — 100 
belonging to the main peak of the local recurrence-time dis- 
tribution function. The transmission parameter is a = 0.17. 
As the events repeat, the stress distribution at the time of 
toppling gets more and more concentrated on the borderline 
value /c = 1. 
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FIG. 10: (Color online) Four types of epicenter sites at the 
center of the figure (blue). Shaded region represents the rup- 
ture zone. 



Next, we examine where in the rupture zone an epicen- 
ter site is located for both cases of the asperity event and 
for the non- asperity event. We classify the site i into four 
types according to the number rii of its nearest-neighbor 
sites which topple during a given avalanche. The cases 
71^=1,2,3 and 4 typically represent an epicenter site lying 
at the tip, corner, boundary and interior of the rupture 
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FIG. 11: (Color online) The rate of the type of epicenter 
sites, for the asperity events (a) and for other non-asperity 
events (b), plotted on a semi- logarithmic scale. The abscissa 
m represents the number of nearest-neighbor sites contained 
in the rupture zone: See Fig. 10. The transmission parameter 
is a = 0.20. 



zone, respectively: See Fig. 10. In Fig.ll, we show on a 
semi-logarithmic plot the rate of each type of an epicen- 
ter site of large events with s > Sc = 100, for the asper- 
ity events (a) and for the non- asperity events (b). The 
transmission parameter is a = 0.20. As can be seen from 
Fig.ll, there is a pronounced tendency that the epicenter 
site is located at smaller n^, particularly at = 1 corre- 
sponding to the tip of the rupture zone. This tendency 
is more enhanced in the asperity events than in other 
events. It turns out that the interior site with n^=4 can 
never be an epicenter of the asperity events while a small 
number of such events (~ 10%) exist for the non-asperity 
event. 

Clear difference exists between the asperity and the 
non- asperity events in other quantities as well. In Fig. 12, 
we show for the case of a = 0.2 the temporal varia- 
tion of the frequency of both the asperity events and 
the non-asperity events observed in a typical run of a 
512 X 512 lattice. The frequency is defined here as the 
number of each type of events occurring in a time bin 
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of one T* (= 1 — 4a = 0.2). As can be seen from the 
figure, the frequency of the non-asperity events fluctu- 
ates around a constant value, whereas that of the as- 
perity events exhibits much wilder temporal variation, 
occasionally vanishing altogether. In other words, there 
appear in turn an active period and a calm period in the 
asperity events, while no such distinction seems to exist 
in the no n- asperity events. 
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FIG. 12: (Color online) The time variation of the frequency of 
the asperity events and of the non-asperity events in a typical 
run of a 512 x 512 lattice. The transmission parameter is 
a = 0.2. The frequency is defined by the number of each type 
of events which occur in a time bin of one T* (= l — Aa — 0.2). 

The asperity events tend to cluster not only in time 
but also in space. In Fig. 13, we show for the case of 
a = 0.2 typical spatial distributions of both the asperity 
events and the non-asperity events on a 256x256 lattice. 
The rupture zone of both types of events, which occur 
in a time period of one T*, are shown together with the 
epicenter site of the asperity events. Fig. (a) corresponds 
to the midst of the "asperity period" where the frequency 
of the asperity events takes its maximum, while Fig.(b) 
corresponds to the midst of the "non-asperity period" 
where the frequency of the asperity events takes its mini- 
mum. As can be seen from Fig. 13(a), the asperity events 
tend to cluster spatially: In particular, the epicenter site 
of the asperity event, which tends to lie at the tip or at 
the corner of the rupture zone, is contiguous either to the 
rupture zone of other asperity events or to the area of no 
(small) events, but not contiguous to the rupture zone of 
non- asperity events. Since the multi-toppling rarely oc- 
curs in the model, this observation means that the stress 
the epicenter site of the asperity event receives from its 
neighboring sites tends to be nearly equal to a. There 
also exists a tendency that the asperity cluster lies in the 
interior of the lattice rather than near the boundary, as 
was already noticed in Ref . [lit • 

In the "non- asperity period", by contrast, large events 
are dominated almost exclusively by the non-asperity 




FIG. 13: (Color online) The spatial distribution of large 
events in a typical run of a 256x256 lattice. The rupture 
zone of both the asperity events and the non-asperity events, 
which occur in a time period of one T*, are shown. Fig. (a) 
corresponds to the midst of the "asperity period" where the 
frequency of the asperity events takes its maximum, while 
Fig.(b) corresponds to the midst of the "non-asperity period" 
where the frequency of the asperity events takes its minimum. 
The transmission parameter is a = 0.2. Dark blue area rep- 
resents either only small events of its size s < Sc = 100 or no 
event occurring during the time interval of T*. Thin orange 
(thick gray) area represents asperity (non-asperity) cluster of 
its size s > Sc = 100. The red bold line represents the rup- 
ture zone of each large event, while green square represents 
an epicenter site of the asperity event. 



events: See Fig. 13(b). Since the number of the non- 
asperity events is more or less kept constant in time as 
can be seen from Fig. 12, the size of the non- asperity clus- 
ter tends to get larger in the non-asperity period, which 
seems consistent with the pattern observed in Fig. 13(b). 

In the following, we shall concentrate on the asperity 
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events. An epicenter site topples at the beginning of a 
given event releasing its stress to zero. If the epicenter 
site i is of the type (n^ = 1, 2, 3, 4), it receives the stress 
from its nearest-neighbor sites times during this asper- 
ity event. Since any toppling in the asperity event occurs 
close to the borderline stress /c = 1, the stress value the 
epicenter site receives during this event is equal to n^a. 
As mentioned, in asperity events, the same site i is likely 
to be an epicenter of the next event again, i.e., reaches 
the threshold stress fc = ^ earlier than any other site in 
the rupture zone of this event. This might sound a bit 
surprising at first, since, as shown in Fig.ll, an epicenter 
site tends to be located at the tip position of the rupture 
zone with = 1, which means that the epicenter site is 
in a low stress state just after the preceding event. How- 
ever, such an epicenter site rapidly develops its stress not 
only by stress loading but also by receiving stress via the 
topplings of its 4 — nearest-neighbor sites not belonging 
to the rupture zone of the preceding event. Here recall 
that the stress transfer the epicenter site i receives from 
the nearest-neighbor sites during the interval time tend 
to be nearly equal to a. 

Then, one can see the reason why the same epicen- 
ter site becomes an epicenter of the next asperity event 
again. The stress the epicenter site i receives from the 
one-site events (or other asperity events) during the in- 
terval time is equal to {4 — ni)a so that the time needed 
for the site i to reach the threshold stress /c = 1 is 
T = 1 - Uia - (4 - ni)a = 1 - 4a = T*, which is just the 
interval time of the asperity events observed in the local 
recurrence-time distribution. 

Next, let us consider the stress at an arbitrary site i 
in the rupture zone other than the epicenter site. We 
suppose that, during a given asperity event, among 
rii nearest-neighbor sites topple before the site i and the 
remaining — nearest-neighbor sites topple after the 
site i. By definition, 1 < < n^. Then, the stress at 
the site i just after the asperity event is {ui — mi)a. The 
stress the site i receives during the interval time from the 
sites not belonging to the rupture zone of the preceding 
event is {A — ni)a. Hence, the time needed for the site i to 
reach the threshold stress /c = 1 is T = 1 — (n^ — m^ja — 
(4 — = 1 — Aa ^ ruia, which is always greater than 
the corresponding time for the epicenter site estimated 
above, 1 — 4a, since > 1. Thus, an epicenter site of 
the preceding asperity event tends to be an epicenter of 
the next asperity event again. 

One sees that the value of does not matter in the 
proof given above. Therefore, our observation in Fig.ll 
that a tip (or a corner) site tends to be an epicenter 
of the avalanche can not be explained solely from the 
above argument. At the moment, we do not know the 
reason why an epicenter site tends to lie at the tip or the 
corner of the rupture zone rather than in its interior in 
the OFC model, though the tendency is quite distinct in 
our numerical simulation. 
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FIG. 14: (Color online) A typical pattern of the asperity- 
sequence interruption, where one big asperity is divided into 
two parts. The transmission parameter is a — 0.2. Fig. (a), 
(b) and (c) correspond to the time t = to (a), t = to + 0.199978 
(b) and t = to + 0.199979 (c). The region surrounded by the 
red bold line represents the rupture zone of the first event. 



The asperity events repeat over many times, but they 
do not last forever. After all, our model is a spatially 
uniform model so that no interior site can be special in 
the long-time limit. Then, we also study the manner how 
the asperity sequence is interrupted. For the asperity se- 
quence interruption, there seem to be two patterns. The 
first pattern is a detachment of the asperity. Namely, one 
big asperity is divided into two parts, each of which col- 
lapse not at the same time, but at different times, though 
still mutually close in time. An example of this pattern 
is given in Fig. 14. The second pattern is an enlargement 
of the asperity. Namely, a given asperity event gets in- 
volved in other larger event becoming a part of it, which 
can occur when the asperity fails to rupture at the regu- 
lar period T* = 1 — 4a due to some reason. An example 
of the latter pattern is given in Fig. 15. 



SUMMARY AND DISCUSSION 

We studied the properties of the OFC model of earth- 
quakes. In this model, characteristic and critical prop- 
erties coexist in an intriguing manner. We computed 
the magnitude distribution, the recurrence-time distri- 
bution, the stress distribution, and investigated various 
asperity characteristics in some detail. Interestingly, the 
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FIG. 15: (Color online) Another typical pattern of the 
asperity-sequence interruption, where a given asperity event 
gets involved in other larger event becoming a part of it. The 
transmission parameter is a = 0.2. Fig. (a), (b) and (c) cor- 
respond to the time t = to (a), t = to + 0.2000 (b) and 
t = to + 0.4248 (c). The region surrounded by the red bold 
line represents the rupture zone of the first event. 



OFC model exhibits not only well-known critical fea- 
tures like the GR law and the Omori law, but also pro- 
nounced characteristic or periodic features, most notably 
exemplified in the occurrence of asperity-like phenomena. 
We have found that the stress dissipation, i.e., the non- 
conservation of the stress, is essential for the asperity-like 
phenomena to occur. We have also shown that a key in- 
gredient in the asperity formation is a self-organization 
of the highly concentrated stress state. Such a stress con- 
centration immediately explains why the interval time of 
the asperity events is equal to 1 — 4a, and why the same 
site becomes an epicenter in the asperity sequence. The 
asperity-like events in the OFC model closely resemble 
those familiar in seismology in the sense that almost 
the same spatial region ruptures repeatedly with a com- 
mon epicenter site and with a common period. 

Indeed, apparent coexistence of critical and charac- 
teristic features is a notable feature of real earthquake 
phenomena as well. The critical features are most typ- 
ically seen in various power-laws observed in statisti- 
cal properties of earthquakes, e.g., the GR law and the 
Omori law as described. Typical examples of character- 
istic features of real earthquakes may be, e.g., the asper- 
ities observed along the subduction zone in northeastern 
Japan, particularly repeating earthquakes off Kamaishi 



22I , [isj . In seismology, the concept of earthquake cy- 
cles has been used in long-term probabilistic earthquake 
forecasts [l|, [H-iel . 

Thus, the OFC model, though an extremely simplified 
model, may capture some of the essential ingredients nec- 
essary to understand apparent coexistence of critical and 
characteristic properties in real earthquakes. In seismol- 
ogy, the origin of the asperity is usually ascribed to possi- 
ble inhomogeneity in the material property of the crust or 
in the external conditions. We wish to stress here, how- 
ever, that in the present OFC model there is no built-in 
inhomogeneity in the model parameters nor in the exter- 
nal conditions, yet the "asperity" is self-generated from 
the spatially uniform evolution-rule and the model pa- 
rameter. 

As mentioned, the asperity in the OFC model is not a 
permanent one. After all, the model is spatially uniform. 
Meanwhile, we have observed that the asperity events 
tend to cluster both in time and in space. In particular, 
the asperity events exhibit "active" and "calm" periods 
in turn. In its active period, the asperity often exists sta- 
bly over many earthquake recurrences. Thus, one needs 
to keep in mind a possibility that, even in the asperity 
formation of real seismicity, inhomogeneity might play 
only a secondary role. 

Of course, real crust is certainly inhomogeneous, and it 
is important to clarify how such inhomogeneity affects, or 
does not affect, the critical and the characteristic prop- 
erties of earthquakes. In this context, it might be in- 
teresting to study the corresponding properties of the 
inhomogeneous OFC model in comparison with those of 
the present homogeneous OFC model. Among others, 
the question whether the asperity-like phenomena and 
the associated characteristic properties would survive the 
build-in inhomogeneity or not is of particular interest. 
The answer to this question seems to be "yes" , although 
some of the properties change from the ones observed 
for the present homogeneous model and its details often 
depend on the specific parameter the model. The proper- 
ties of such inhomogeneous OFC model will be reported 
in our forthcoming publication [27]. 

The authors are thankful to Prof. N. Kato for useful 
discussion. This study was supported by Grant-in-Aid 
for Scientific Research 19052006. We thank ISSP, Tokyo 
University for providing us with the CPU time. 

APPENDIX 

In this appendix, we show for the asperity event that 
the stress distribution at the time of toppling tends to 
concentrate more on the threshold value /c = 1 as the 
asperity events repeat. 

First, we look at the site j which is the nearest neighbor 
of an epicenter site 0. Let the stress of the site j be 
fj = 1 — a -\- X. We assume < x << 1, which is 
equivalent to assuming that this site topples nearly at its 
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threshold value /c = 1. Among the four nearest-neighbor 
sites of the site j, one is an epicenter 0. The epicenter 
site always topples at / = /c = 1 and distribute the 
stress a to its nearest-neighbor sites including j. Let the 
stress of other three nearest-neighbor sites of the site j 
he fi = 1 - a ^ X ^ SiX {i = 1,2, 3): See Fig.l6. After 
the toppling of the epicenter site 0, the site j reaches 
the stress 1 -\- x > 1 and topples. The three nearest- 
neighbor sites of the site j then reach the stress fi-\-afj = 
1 -\- X -\- ax -\- diX {i = 1, 2, 3) and topple. The site j then 
receives the stress from the sites 1, 2, and 3, reaching 
the stress 3a(l + x + ax) + a{Si -\- 62 -\- Ss)x. As long 
as the site j does not topple again in this asperity event, 
which is usually the case, this is the stress value at the 
site j immediately after this asperity event. After the 
interseismic period T* = 1 — 4a, the stress at the site j 
becomes 1 — a-\-3a{l -\-a-\-5)x with S = {Si +^2 +^3)/3. 
Then, the condition of the stress concentration occurring 
at the site j is 3a{l -\- a -\- 6) < 1. Namely, as the asperity 
events repeat, the stress at the time of toppling gradually 
approaches the threshold value /c = 1 if the condition 
^ < (5c = ^ - (1 + a) is fulfilled. 
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FIG. 16: Illustration of neighboring toppling sites in the 
rupture zone, with its initial stress values. 

Now, we extend the above discussion to more general 
situation where, among the four nearest-neighbor sites 
of the site j, (i) one site topples before the site j and 
three other sites 1,2 and 3 topple after the site j. Here 
we consider the site j to be an interior site with all its 
nearest-neighbor sites contained in the rupture zone of 
the asperity event. Let the stress of the site j be fj = 
1 — a -\- X as above, the stress of the site toppling before 
j be /o = 1 + X + Sqx, and the ones after j be fi = 
1 — a -\- X -\- 5iX {i = 1, 2, 3). After the toppling of the 
site 0, the site j reaches the stress 1 + x + a{l -\- So)x 
and topples. The three neighboring sites 1, 2 and 3 then 
reach the stress fi + afj = 1 -\- x -\- ax -\- 5iX -\- a'^{l-\- 5o)x 
{i = 1,2,3) and topple. The site j then receive the stress 
from the sites 1,2, and 3, reaching the stress 3a(l + x + 
ax^a^x)^3a5x^?>a^5QX. As long as the site j does not 
topple again in this asperity event, which is usually the 



case, this is the stress value at the site j just after this 
asperity event. After the period T* = 1 — 4a, the stress 
at the site j becomes 1 — a + 3a(l + a + a^ + (5 + a^5Q)x. 
Then, the condition of the stress concentration occurring 
at the site j is ^ < = — {1 + <^ + <^^(1 + ^')} where 
5' = So. 
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FIG. 17: (Color online) The threshold value for the stress 
concentration, 6c, plotted versus the transmission parameter 
a in the case of 6' = 0.1. The three processes (i) ^ (iii) are 
described in the text. One can see that the stress concentra- 
tion process, expected at ^ < ^c, is readily operative in most 
cases, since Sc is large except for the process (i) close to the 
conservation limit a — 0.25. 



One can repeat similar discussion to the other cases 
where, among the four nearest-neighbor sites of the site 
j, (ii) two sites and 1 topple before the site j and the 
two other sites 2 and 3 topple after the site j, (iii) three 
sites 0, 1 and 2 topple before the site j and one other site 3 
topples after the site j, and (iv) all four sites topple after 
the site j. In the case (ii), the condition of the stress 
concentration is S < Sc = ^ — {l + a + 2a^(l -\-S')} with 
S = {S2^Ss)/2 and S' = {So^Si)/2. In the case (iii), the 
condition of the stress concentration is S < Sc = ^ — {I -\- 
a + 3a2(l + S')} with S = S3 and S' = {So + (5i + S2)/3. 
The case (iv) is exceptional, where the stress of the site 
j after the event is always zero. The a-dependence of Sc 
is shown in Fig. 17 in the above cases of (i) ~ (iii) with 
S^ = 0.1. In fact, the (5c- value is rather insensitive to the 
(5^- value, since, in the expression of Sc, the coefficient of 
the S'-teim is smaller than the leading term by factor of 
O(a^) with a < 0.25. 

The above analysis indicates that, once each site starts 
to topple more or less at similar stress values close to 
the threshold value /c = 1, this tendency is more and 
more evolved as the asperity events repeat, i.e., the stress 
concentration tends to be self-organized. 
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